#pragma once
#include "elisa/elisa.h"
#include <iostream>
#include <sstream>
#include <algorithm>
#include <string>
#include <vector>
#include <gaol/gaol>

#include "Modeler.h"

#ifdef ELISA_PERFORMER
#include "elisa_chrono.h"
#endif

using namespace std;
using namespace elisa;

/* Problem: the well-known transistor problem described for instance in

  @book{Numerica,
  author    = {P. {Van Hentenryck} and L. Michel and Y. Deville},
  title     = {{N}umerica: a {M}odeling {L}anguage for {G}lobal {O}ptimization},
  publisher = {MIT Press},
  year      = {1997}
  }
*/
class Transistor : public Modeler
{
private:
	Model _model;

	RealConstant G_1_1, G_1_2, G_1_3, G_1_4;
	RealConstant G_2_1, G_2_2, G_2_3, G_2_4;
	RealConstant G_3_1, G_3_2, G_3_3, G_3_4;
	RealConstant G_4_1, G_4_2, G_4_3, G_4_4;
	RealConstant G_5_1, G_5_2, G_5_3, G_5_4;
	RealVar X[9];
	RealVar Y[3];
	RealConstraint c1, c2, c3, c4, c5, c6, c7, c8, c9, c10, c11, c12;

public:
	Transistor()
		: G_1_1(0.485), G_2_1(0.369), G_3_1(5.2095), G_4_1(23.3037), G_5_1(28.5132),
		  G_1_2(0.752), G_2_2(1.254), G_3_2(10.0677), G_4_2(101.779), G_5_2(111.8467),
		  G_1_3(0.869), G_2_3(0.703), G_3_3(22.9274), G_4_3(111.461), G_5_3(134.3884),
		  G_1_4(0.982), G_2_4(1.455), G_3_4(20.2153), G_4_4(191.267), G_5_4(211.4823),
		  c1(X[0]*X[2] == X[1]*X[3]),
		  c2(X[0]*X[1] == 1-Y[0]),
		  c3(Y[1]      == Y[0]*X[2]),
		  c4(Y[2]      == Y[0]*X[3]),
		  c5(Y[1]*(exp(X[4]*(G_1_1 - G_3_1*X[6]*0.001 - G_5_1*X[7]*0.001)) - 1) == G_5_1 - G_4_1*X[1]),
		  c6(Y[1]*(exp(X[4]*(G_1_2 - G_3_2*X[6]*0.001 - G_5_2*X[7]*0.001)) - 1) == G_5_2 - G_4_2*X[1]),
		  c7(Y[1]*(exp(X[4]*(G_1_3 - G_3_3*X[6]*0.001 - G_5_3*X[7]*0.001)) - 1) == G_5_3 - G_4_3*X[1]),
		  c8(Y[1]*(exp(X[4]*(G_1_4 - G_3_4*X[6]*0.001 - G_5_4*X[7]*0.001)) - 1) == G_5_4 - G_4_4*X[1]),
		  c9(Y[2]*(exp(X[5]*(G_1_1 - G_2_1 - G_3_1*X[6]*0.001 + G_4_1*X[8]*0.001)) - 1) == G_5_1*X[0] - G_4_1),
		  c10(Y[2]*(exp(X[5]*(G_1_2 - G_2_2 - G_3_2*X[6]*0.001 + G_4_2*X[8]*0.001)) - 1) == G_5_2*X[0] - G_4_2),
		  c11(Y[2]*(exp(X[5]*(G_1_3 - G_2_3 - G_3_3*X[6]*0.001 + G_4_3*X[8]*0.001)) - 1) == G_5_3*X[0] - G_4_3),
		  c12(Y[2]*(exp(X[5]*(G_1_4 - G_2_4 - G_3_4*X[6]*0.001 + G_4_4*X[8]*0.001)) - 1) == G_5_4*X[0] - G_4_4),
		  _model()
	{
		for(int i=0; i<9; i++) {
			string name("X["), e("]");
			ostringstream oss;
			oss << (i+1);
			name += oss.str();
			name += e;
			new( X + i ) RealVar(Interval(0,10), name);
		}
		for(int i=0; i<3; i++) {
			string name("Y["), e("]");
			ostringstream oss;
			oss << i;
			name += oss.str();
			name += e;
			new( Y + i ) RealVar(Interval::Universe, name);
		}
	}

	elisa::Model & model()
	{
		_model = c1 && c2 && c3 && c4 && c5 && c6 && c7 && c8 && c9 && c10 && c11 && c12;
		return _model;
	}

	void print()
	{
		for(int i=0; i<9; i++) {
			cout << '\t' << X[i].getName() << " = " << X[i].domain() << endl;
		}
		for(int i=0; i<3; i++) {
			cout << '\t' << Y[i].getName() << " = " << Y[i].domain() << endl;
		}
	}

	RealVar & varX()
	{
		return X[0];
	}
	RealVar & varY()
	{
		return Y[0];
	}
};
